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Abstract 

Using a Monte Carlo coarse-graining technique introduced by Binder et al, 
we have explicitly constructed the continuum field theory for the zero- 
temperature triangular Ising antiferromagnet. We verify the conjecture that 
this is a gaussian theory of the height variable in the interface representation 
of the spin model. We also measure the height-height correlation function and 
deduce the stiffness constant. In addition, we investigate the nature of defect- 
defect interactions at finite temperatures, and find that the two-dimensional 
Coulomb gas scenario applies at low temperatures. 

PACS numbers: 05.70.Jk, 05.50.+q, 05.10.Ln, 02.70.Lq 
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I. INTRODUCTION 



In recent years, there has been considerable interest in the study of classical spin systems 
with critical ground states. One of the best studied of these is the triangular-lattice Ising 
antiferromagnet (TIAFM) |T[. The critical behavior of these models can be understood 
based on an interface representation and an "effective" field theory which is gaussian in the 
height variable . In this paper, we present an explicit construction of the effective field 
theory based on the study of Monte Carlo cell distribution functions [Q of the TIAFM. To 
accomplish this, a height mapping is applied to the system and this variable is then coarse 
grained to obtain a continuum field theory. Our interest in this model was stimulated by the 
observation of anomalously slow dynamics in the compressible TIAFM ||. In that model, 
the supercooled state exhibits an ergodicity-breaking transition which is reminiscent of the 
structural-glass transition. This glassy behavior is believed to arise from the interaction 
between strings and vortices 0, the topological defects present in these models at finite 
temperature. In the pure TIAFM model, the gaussian theory implies that the vortices 
interact as charges of a 2-D Coulomb gas || . We have analyzed the defect-defect correlation 
function at finite temperatures to investigate the nature of the defect interactions. The 2-D 
simulations are consistent with the Coulomb gas picture at low densities of defects. 

The motivation behind this numerical study was (1) to test the Monte Carlo cell dis- 
tribution function technique in a model where the effective field theory is well established 
and (2) to establish a framework for the construction of effective field theories and effective 
defect-defect interactions in models such as the compressible TIAFM where no adequate 
field theory description exists. 

The paper is organized as follows. In section II we describe the cell distribution function 
technique. In section III we present results from the coarse grained free-energy functional 
and compare our results to the gaussian theory. In section IV we verify the gaussian theory 
from a study of the height-height correlation function and in section V we describe the study 
of the defect-defect correlation functions. Section VI presents our conclusions and directions 
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for future work. 



II. CELL DISTRIBUTION FUNCTIONS AND COARSE GRAINED 
FREE-ENERGY FUNCTION ALS 

Binder and coworkers P,|T0|| have studied cell distribution functions of the three di- 
mensional Ising model on a cubic lattice and constructed coarse-grained Ginzburg-Landau 
Hamiltonians. In this section, we review this technique and present the results of its appli- 
cation to the triangular Ising ferromagnet as an example. We will describe the application 
of this technique to the zero-temperature TIAFM in the next section. 

Given a microscopic Hamiltonian such as the Ising model, Monte Carlo(MC) methods 
can be used to sample the distribution functions Pi({sj}) of the coarse grained variables 

where Si are the original microscopic spin variables, such as the Ising spins on the original 
lattice, and L is the cell size. The distributions Pi({sj}) are assumed to be of the form e~^ T 
with T having the Ginzburg Landau(GL) form in terms of the coarse-grained variables {sj}. 
This assumption is expected to be valid when L » a(lattice spacing) but much smaller 
than the correlation length, such that the coarse grained variables do not fluctuate rapidly 
from cell to cell. If one could sample the total distribution function in MC simulations, then 
this connection could be exploited for the explicit construction of the GL Hamiltonian by 
simulating the microscopic model. Sampling the total distribution function is essentially an 
impossible task and therefore, we have followed Binder jT(|] in studying the two simplest 



reduced distribution functions, the single cell and the joint, nearest-neighbor, two-cell dis- 
tribution functions which are then parameterized by the GL form. To illustrate how the 
cell distribution function method works, we briefly describe its application to the triangular 
Ising ferromagnet in zero magnetic field, where we choose the model system size to be 80x80 
with periodic boundary conditions and the cell sizes L to be 4, 8, 10, 16, 20 and 40. 
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The microscopic Hamiltonian is: 



Wising — —J E ^iSi 



<w> 

The total distribution function P L ({.Si}) is assumed to be: 

PUisi}) = |e-^«^» (1) 
where Z is the partition function and the GL Hamiltonian has the form 

^gl({^}) = + u Lst) + £ - s 3 f . (2) 

i <ij> 

The two reduced distribution functions, which are amenable to numerical calculations, are 
the two-cell joint distribution function P L (s i ,s J -) and the single-cell distribution function 

P L (s l ,s j )= f n d Sl P L ({ Sl }) (3) 

<ij> 

^L(si) = J P L {s hSj )ds 3 (4) 

The parameterization of the effective Hamiltonian is accomplished by sampling Pt(sj, Sj) 
and Pl(sj) 

s i) = Y< e M~Msi " s if + Vb(si) + Vi( Sj ))} (5) 

ftW = (6) 

where Z' and Z" are normalization factors, and 

V L (si) = -r L s1 + u L s\ 

The rL and ul obtained from Pl(si, Sj) and Px,(sj) respectively, can be different since 
the effect of the gradient term has been integrated over in the single-cell distribution func- 
tion. Similar considerations imply that these coefficients can be different from rL, ul and 



cl in Eq. 0. For small values of c~l, the coupling between different sites is small and all 
of these coefficients are expected to be approximately equal [TIJ. The coefficients r L , ul 
and cl, as functions of temperature (T) and cell size (L), can be estimated by fitting the 
measured distribution function PL{si, Sj) to the model form (Eq. It is well known that a 
temperature-driven second order transition exists for this Ising model and the exact result 
for T c (oo) = 3.6 PHI . Fig. p] shows as a function of temperature for different cell sizes. 
There is clear evidence for a temperature-driven second order transition at T C (L) where the 
sign of the ratio changes. The L dependence of T C (L) agrees with finite size scaling pre- 
dictions H|]. The variation of with cen size L (Fig. ||) shows temperature-dependent 
behavior as the size of the cells approach infinity. For temperatures above T c (oo), the ratio 
approaches a non-zero fixed point as L approaches infinity. For all other temperatures, this 
ratio approaches zero. The sampling of the single-cell distribution function Pi{si) (Eq. ^) 
leads to similar results for tl and ul, which is consistent with the small value obtained for 
cl- These results support the assumptions that went into the construction of J-'glUsi})- 

The results of this section will be contrasted with the frustrated antiferromagnetic case 
in the next section. The application of the distribution-function technique to the TIAFM 
is based on the mapping to an interface model and height variables. The coarse graining of 
these height variables and the GL parameterization are described in the next section. 



III. RESULTS OF COARSE GRAINING AND COMPARISON WITH THE 

GAUSSIAN THEORY 

A. Height Variables and coarse graining at T=0 

One interesting property of the TIAFM is that the ground state ensemble has a one-to- 
one mapping onto an interface model @J^,§|. This provides a simple way of studying the 
properties of zero-temperature TIAFM. 

The discrete variables in the interface model are height variables {zi}, which are defined 
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on each site of the triangular lattice. The mapping from the TIAFM to an interface model 
is realized by the mapping from spin variables {Si} to height variables {zi}. Specifying the 
height variable to be zero at a chosen site, the following rules (Fig. |3[), provide a unique 
mapping from {Si} to {zi}: 

(1) Along dl direction, Az = — 1 for opposite spins and Az = +2 for same spins. 

(2) Along d2 direction, Az = — 1 for opposite spins and Az = +2 for same spins. 

(3) Along d3 direction, Az = +1 for opposite spins and Az = —2 for same spins. 

These height assignments are unique up to the choice of origin as long as we restrict ourselves 
to Ising configurations in the ground state ensemble, which is the set of states that do not 
have any completely frustrated plaquettes (all spins same). The microscopic Hamiltonian 
in term of these height variabls {z{{ can be written directly from the mapping as [0 

H{z t } = -J E ( 3 -21^-^1) 

<ij> 

The symmetries of the original Hamiltonian are reflected in 7i{z{}. These include Ising up- 
down symmetry, and sub- lattice symmetry. In addition, 7i{zi} has a symmetry with respect 
to global, discrete shifts of the height variable. A detailed discussion of these symmetries 
can be found in Blote et al. ||. An interesting property of this mapping is that for any 
ground state the height variable modulo 3 is the same for all spins on the same sub-lattice. 

Starting from {z/}, we can define a new set of height variables {H^} which are situated 
on the dual lattice, by averaging the three neighboring z\ on each triangular plaquette: 
H/j, = \{z\ + z 2 + z 3 ) where z 1: z 2 , and z 3 are on the three vertices of the fj, th plaquette [12| 



as shown in Fig. |3|. Since the sum of z\, z 2 , and z 3 are multiples of 3, the {H^} is a set of 
integers. 

We apply the coarse graining technique to the height variables {H^}. The coarse-graining 
cell is chosen to be a rhombus with the linear size L and the coarse-grained height variables 
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{hi} are defined as the average of all the within the cell. 

hi = ± £ H, 



:th 



where N is the total number of variables within the i cell. The coarse-graining cell 
is depicted schematically in Fig. || As L approaches infinity, the {hi} become continuous 
variables. The effective field theory is based on these continuous variables {hi}, which 
describe the roughness of the interface. In the rough phase of the interface model, the 
average tilt of the interface, which is defined as the difference between < h m > and < h„ > 
with m th cell and n th cell separated by the system size, is zero. For an arbitrary spin 
configuration, a tilt can be frozen in ||. We will avoid such configurations in our simulations. 
Also {< hi >} may be non-zero for a general choice of height origin and the Hamiltonian 
should only depend on the fluctuation h— < hi > 0. 

We performed Monte Carlo simulations on a system of size 600x600 shaped as a rhombus 
and looked at cell sizes 10, 15, 20, 30, 40, 50, 60, 75, 100. We set the initial configuration to 
be one of the flat states fll2| , where {H^} is uniform in space. The dynamics we use in the 



simulation is single spin-flip Metropolis algorithm at zero temperature, where the energy is 
not allowed to increase and thus not all spins are flippable. The first 20,000 configurations 
are discarded to ensure that measurements are taken in equilibrium. The sampling is done 
every ten MC steps and for a total of 6,000 configurations. Since we use single spin-flip MC 
dynamics, the average tilt is unchanged from its initial value of zero during the MC runs 
|T3|j . In the following sections, we will use hi to represent hi— < hi >. 



B. Results and comparison with the gaussian theory 

It has been conjectured |H§ that in terms of the height variables, the continuum theory 
for the zero temperature TIAFM is gaussian: 

F{h{r)) = J dr C - (V/i(f)) 2 (7) 
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with a value of c = 7r/9. This stiffness constant was identified by Blote et al through 
the correspondence between the exact calculation of the spin-spin correlation function in 
the spin model |T^] and the discrete height-height correlation function in real space in the 



interface model |J . The stiffness constant was also obtained numerically by Zeng and Henley 
from measurements of the discrete height correlation function in Fourier space [E2|. In this 



section, we deduce the stiffness constant by an explicit construction of the coarse-grained 
free energy. 

Assuming a Ginzburg-Landau form for the Hamiltonian in terms of the continuous height 
variables, 

F(h(r)) =J dr{\ W? + ~ (V%1) 2 } (8) 

we proceed to obtain a and c from our simulations by parameterizing the joint distribution 
function of the coarse-grained height variables hj). 

As in the triangular ferromagnet, we study the nearest-neighbor joint distribution func- 
tion of heights, 

P L (h t ,h J )= f J] dhP^ih}) (9) 

<ij> 

where Pi({/i;}) is the total cell distribution function. We parameterize the joint distribution 
function as 

Pl{K hj) = ^ ex V {-[c L (h t - htf + a L {h\ + h))}} (10) 

where Z' is a normalization factor. 

Anticipating a critical point at T=0, where cll — > 0, we have neglected the quartic term 
in Eq. |10[ This will be justified a posteriori from our numerical study. From the joint 
distribution function, we can extract and cl as in the ferromagnetic Ising model. 

We examine a typical Pi,(h,h*) for cell size L=50. One-dimensional cuts of Pi(hi,hj) 
at different fixed values hj = h* are shown in Fig. f|. The peak and shape of the cuts can 
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be deduced from the distribution function (Eq. [10]) which is a product of two Gaussians, 
exp(— aih 2 ) and exp(— cl{H — h*) 2 ). When h* differs from by a large amount, the plots in 
Fig. |] show that the actual peak of the cut is much closer to h* than to 0, which implies that 
the exp(— cl(/i — h*) 2 ) is the dominant term in the product of Gaussians and the exp(— aLh 2 ) 
term acts more as a pre-factor modulating the amplitude of the peak. The width of the cut 
is also seen to be determined mainly by cl- Fig. |] indicates that is small compared 
with cl- We also analyze the importance of higher order gradient terms by measuring the 
non-gaussian parameter (gi) in the distribution of v = \{hi — hj): 

< V* > L 

9l = 1- o ~ 2 o (11) 

3 < v 2 > L 

We find gi to be much smaller than 1, which is consistent with the assumed form of the 
gradient term in the joint distribution function. Extracting quantitative information about 
the coefficients at, and cl from fits to the distribution functions, turns out to be difficult 
because of the essential two-dimensional nature of the distribution function Pi,{hi, hj)(cL is 
large). Instead, we resort to measurement of moments for extracting quantitative informa- 
tion. 

The parameters a^ and cl can be related to the moments of {hi}. We have calculated 
various moments of these coarse-grained variables directly from the MC simulations. In 
general, the parameterization of the joint distribution function (such as Eq. [H]) provides 
a connection between these moments and the coefficients of the GL Hamiltonian. For the 
model of Eq. [l(], the relation between ai and Cl and the moments can be shown to be: 

1 

aL ~ 2(< h 2 > + < hji2 >) 

= < hih 2 > ( s 

° L 2(< h 2 > + < hih 2 >)(< h 2 >-< h x h 2 >) 1 ' 

implying that 

cl < hih 2 > ^ 



a L (< h 2 > - < hih 2 >) 



Here hi and h,2 refer to coarse-grained height variables on nearest-neighbor cells separated 
by the coarse-graining cell size L. These relations, between the parameters and cl and the 
moments, also provide the connection between cll and cl and the corresponding parameters 
a and c in continuum field theories such as in Eq. |8]. Using the known expressions for the 
correlations functions in the continuum model in two dimensions 111511, these relations are: 



a = 4a L L- 2 (15) 

and 

e = 2u> !^Mc i (16) 

7T 

The factor wq = Vo/2 arises from the volume per unit cell in the triangular lattice. The 
prediction from the gaussian theory is a = and c = tt/9. We expect that, if the gaussian 
theory is correct, then a should approach zero and c should remain constant as L— >oc. 
Because of the connection to the continuum theory(Eq. [L5|), it is natural to fit a^/L 2 , to a 
second-order polynomial in 1/L. The results of the fitting are shown in Fig. [5| The value 
of a can be deduced from the constant term extracted from the fitting and is found to be 
6.183* 10~ 6 implying that the quadratic term in the Ginzburg-Landau Hamiltonian becomes 
negligible as L— >oc. It is difficult to extract the value of c from the measured values of cl 
because of the fact that the errors increase with L | 16| . If, instead, we analyze cl/L 2 , the 
error bars actually decrease with 1/L as seen from Fig |6], and this provides a more precise 
way of determining the value of c. As shown in Fig. || the data for Cl/L 2 can be fit very well 
to a second-order polynomial in 1/L with vanishingly small coefficients of the constant and 
linear term in 1/L. The value of cl obtained from this fit is 0.622. The standard deviations, 
calculated from the appropriate higher moments of hi and hj based on Eq. [1^ and Eq. [L3|, 
are plotted as error bars. 

These results for the parameters obtained from our analysis of the moments are com- 
pletely consistent with the qualitative behavior deduced from the joint distribution functions, 
Pi(hi, hj), shown in Fig. |. 
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The results shown in Fig. |5| and Fig. || verify the gaussian nature of the continuum theory. 
The numerical value of c extracted by applying Eq. [16] to the fitted value of Cl is 0.392, which 
is in very good agreement with the conjectured stiffness constant of 7r/9 = 0.349 [@,[3]]. 

IV. HEIGHT-HEIGHT CORRELATION FUNCTION 

In this section we will briefly present our calculation of the correlation function of height 
variables in Fourier space, i.e., the power spectrum of the height variable, which can be 
used to deduce the stiffness constant as in the investigation of Zeng and Henley |12| . In the 



next section we further extend this technique to study the defect-defect interactions at finite 
temperatures by measuring the defect density correlation function. 

As implied by Eq. [7|, the stiffness constant c of the interface model can be directly related 
to the correlation function of height variables {hi} in Fourier space. 

T{h{r)} = Jdf°- (Vh(r)f = £ \ |«f \h(q) | 2 (17) 

q 

Thus, the correlation function < \h(q)\ 2 > is: 

<I%T>=-^ (18) 
c \q\ z 

The stiffness constant c can be extracted from a knowledge of < \h(q)\ 2 >. Numerical 
measurements of < \h(q)\ 2 > are, however, arduous, and instead, we study the correlation 
function of the discrete height variables, < \H(q)\ 2 >, which should resemble < \h(q)\ 2 > 
for small q except for a factor of y/3/2 coming from the volume per lattice site. Therefore, 
the < \H(q)\ 2 > for small q can be written as: 

<\H(Q)\ 2 >=^ (19) 
c \q\ z 



where 

2 



Vs' (20) 
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We simulate a system of size 240x240. We start with the same flat states as the ones 
used in the coarse-graining simulations. After discarding the first 20,000 configurations, the 
simulations are sampled at every 10 MC steps for a total of 8,000 configurations. For each 
sampled configuration we calculate \H(q)\, the two-dimensional Fourier transform of the 
discrete height variables {H^} (cf section IIIA) and average \H(q)\ 2 over the sampled MC 
configurations. 

For simplicity, one-dimensonal cuts of < \H(q)\ 2 > _1 along different directions in q-space 
are shown in Fig. [7[ As seen from Fig. 0, different cuts collapse on top of one another only 
for small values of q 2 . This is due to the fact that the correlation function is isotropic for 
small values of q because of the six-fold rotational symmetry of the triangular lattice. As 
seen from the insert, the anisotropy starts to be significant when q 2 > 1.6. We have fitted 
< \H(q)\ 2 > _1 to a second-order polynomial in q 2 , restricting the fitting region to q 2 < 0.2. 
A value of c = 0.414 ± 0.016 was extracted from the linear term of the fit. The non-linear 
terms were small, and an estimate of the error on c was obtained from the standard deviation 
of the linear term extracted from fits to different segments of the curves (q 2 < 0.2). From 
Eq. the stiffness constant can be deduced to be c = 0.360, compared with 0.392 from 
our previous real space coarse-graining approach and | = 0.349 from Blote et al 0. 

The similarity of the results obtained from the real-space analysis and the power- 
spectrum analysis is quite remarkable and underlines the strength of these numerical tech- 
niques. A comparison of the two different analysis techniques, however, makes it clear that 
the power-spectrum analysis is simpler to implement and is the more desirable technique. 
The question of how well either of these techniques will fare when applied to a model with 
general non-linear interactions, where very little is known analytically, is still open. 

V. DEFECT INTERACTIONS 

At finite temperatures, there are completely frustrated plaquettes which can be viewed 
as vortex excitations in the height variable |fJQ] . There are two types of vortices correspond- 
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ing to the two different orientations of the frustrated triangular plaquettes. These vortex 
excitations interact like charges via a two-dimensional Coulomb potential ||||. 

An alternative to the gaussian description is, therefore, a Coulomb gas description in- 
volving electric and magnetic charges corresponding to spin waves and vortices, respectively 
This connection suggests that extracting effective defect-defect interactions could be 
an alternative to constructing coarse-grained free energies. In complicated models and espe- 
cially in analyzing dynamics, this might be the more viable alternative. We therefore wanted 
to numerically extract the effective vortex-vortex interaction in the TIAFM and verify the 
Coulomb gas scenario. 

In the regime of low defect-density, where mean field theory is expected to hold, the effec- 
tive interactions can be related to the density- density correlation function [fTTfl . The appro- 
priate mean-field theory for our system with positive and negative charges is Debye-Hiickle 
theory |17] which predicts that the charge-density correlation function has the following 
form: 



<|p(g)| 2 >=f -^— 2 (21) 

47T q z + KT 



where 

k 2 = 4tt(< n + > + < n_ >) (22) 

is the square of an inverse Debye-Hiickle screening length. < n± > are the average densi- 
ties of positive and negative charges and p{q)=n + (q) + n_(q), where n±(q) are the Fourier 
transform of the real space charge density functions. 

If the vortices in our model behave as a Coulomb gas, then we expect to find this behavior 
of the charge- density correlation function at low temperatures where the vortex density is 
low. 

We studied the charge density correlation function in a manner similar to that em- 
ployed in the previous section for studying the height-height correlation function. The MC 
simulations were performed on a 120x120 system and after discarding the first 10000 config- 
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urations, samples were taken every 5 configurations. The correlation function, < |p(<?)| 2 >, 
was measured at several temperatures. 

The simulation results are shown in Fig. |8| along with the fits to the Debye-Hiickle 
form(Eq. ^1]). The figure shows a one- dimensional cut of < |p(g)| 2 > along a particular 
direction. In contrast to the height power spectrum, the defect density power spectrum 
was found to be isotropic within the range of q studied. The Debye-Hiickle theory models 
the data well for all temperatures except the lowest temperature, T = 0.6, which is shown 
in the insert. The lowest temperature data is suspect because the defect density is low 
and statistics are difficult to obtain. The isotropy of the < |p(g)| 2 > data imply that the 
screening length and, therefore, the interactions between the defects on the triangular lattice 
is isotropic down to length scales of the order of the lattice spacing. 

From the fits to the Debye-Hiickle theory the inverse square of the screening length, 
K, 2 can be extracted and compared to the analytic form given in Eq. (^) which uses as 



an input the measured defect density. In Fig. |9|, the fitted values, k 2 and the analytic 
values, k 2 . are plotted versus temperature. Except at the highest temperature, T=2.0, the 
agreement between the two different estimates is very good. Interestingly, the temperature 
at which the two k's differ significantly from each other is also where the defect density 
differs from a simple Arrhenius prediction corresponding to an activation energy of 4J, 
the single-defect creation energy. The observation of the large defect density at T=2.0, 
< n + > + < n_ >= 0.110, indicates that the defects are covering the whole lattice at this 
temperature and the single-defect picture becomes inappropriate. 

These numerical studies show that it is possible to extract effective defect-defect in- 
teractions from measurements of the defect-density power spectrum and verifies that the 
interaction, in the TIAFM, is Coulombic. A numerical measurement of the force between 



defects has been made via a study of the dynamics of the defects in TIAFM [[L8|] . This study 
could not confirm the existence of a Coulomb force from the behavior of the defect density 
at long times, however, it was argued that the force could not be falling of more slowly than 
1/r. 
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VI. CONCLUSION 



Based on a mapping to an interface model, we have constructed a coarse-grained Hamil- 
tonian for the zero-temperature TIAFM by studying the joint distribution function of the 
coarse-grained height variables. Our numerical work confirms that the effective field theory 
is gaussian in terms of the coarse-grained height variables and we obtain a numerical esti- 
mate of the stiffness constant which is in a good agreement with the analytical prediction. 
The connection between the parameters in joint distribution function and the moments of 
height variables eliminates the arduous work of directly relating the parameters of joint 
distribution function to those of the total distribution function. We also directly measured 
the power spectrum of the height and obtained a value of the stiffness constant, which is 
very close to the analytic prediction and the estimate from the real-space coarse-graining. 
This suggests that measurement of power-spectrum might provide a simple route towards 
the construction of effective field theories. We also demonstrated that effective defect-defect 
interactions can be extracted from numerical studies of the defect-density power spectrum.. 

The compressible TIAFM |19|j20[| is currently being studied using these techniques. The 
effective field theory and defect-defect interactions in this model are of interest for the study 
of alloys with the elastic interactions, and for the study of glassy dynamics ||. 
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FIG. 2. The ratio CL 1^2 , calculated from the joint distribution function for different tempera- 



tures, plotted as a function of L . The curves can be divided into two groups. For those with 



temperatures higher than T c (oo)=3.6 (identified in the figure), the curves converge to a finite value. 



For those with temperatures lower than T c (oo), the curves converge to zero. 
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FIG. 3. (a) The three nearest-neighbor directions dl, d2, d3, and the 6 nearest neighbors (black 
dots) of discrete height variables Zj(open circle) mapped onto the rhombus-shaped lattice used in 
the simulation, (b) The coarse graining cell of size L = 2. The vertices (black dots) of lattice are 
{z^ and the centers of the triangular plaquettes (open square) are {H^}, each of which is defined as 
the average of the three surrounding z% variables. The coarse-grained height variable h\ is defined 
as the average of all the variables within this cell. 
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FIG. 4. The cuts of h*) for cell size L=50 and different values of h*, demonstrating 

dominance of the gradient term(c/. text) . 
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FIG. 5. Plots of cll = a^/L 2 for different cell sizes. The open circles are the data points and 
the solid line is the second degree polynomial fit. The error bars are obtained from the standard 
deviations calculated from appropriate higher moments of hi and hj . The fitting of / L 2 leads to 
a constant term with a value of 6.183E-6, which implies that the quadratic term in the field theory 
a (cf. Eq. |l5| ) becomes vanishingly small as L — »oo. 
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FIG. 6. Plots of cl/L 2 for different cell sizes. The open circles are the data points and the solid 
line is the second degree polynomial fit. The error bars are the standard deviations calculated from 
appropriate higher moments of hi and hj. The fitting of cl/L 2 leads to a quadratic term with a 
coefficient 0.622, and vanishingly small constant and linear terms, which implies that the gradient 



term of the continuous theory(cf. Eq. 16), c approaches 0.392 as L— >oo. 



23 



40 




0.0 0.4 0.8 1.2 1.6 ,6 s 




2 4 6 8 



FIG. 7. Plots of < \H(q)\ 2 > _1 vs. \q) 2 . Different data sets are from the cuts of two-dimensonal 
data along different directions in q-space. The correlation function shows significant anisotropy for 
q 2 > 1.6. The value of c = 0.414 ±0.016 is obtained from the fitting of curves in the region q 2 < 0.2 
to a second-degree polynomial in q 2 . The stiffness constant, c, deduced from this fit is 0.360. 
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FIG. 8. Plots of the defect density correlation function at T=1.0, 1.2, 1.5, 2.0. The data points 
are taken from cuts of two-dimensional data along one chosen direction since we did not observe 
any anisotropy (c/text), and the solid lines are from the fittings to a Debye-Hiickle form(Eq. |2l| ). 
The fittings lead to at different temperatures. The Kj and the analytical values of k% extracted 
from the total defect density for the same temperature (cf. Eq. |22] ) are shown in Fig. |9|. The total 
defect density is 0.030, 0.051, 0.085, 0.110, respectively, for T=1.0, 1.2, 1.5, 2.0. The inset shows 
the fitting at the lowest temperature, T=0.6 where the defect density is 0.003 and the results suffer 
from a lack of statistics 
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FIG. 9. Plots of K 2 f , the inverse screening length obtained from fitting compared to k 2 , the 



value obtained from Eq. 22. The two curves deviate significantly only at the highest temperature 



where the defect picture itself starts to break down (c/text). 
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